A microscopic quantum dynamics approach to the dilute condensed Bose gas 



Thorstcn Kohlcr and Keith Burnett 
Clarendon Laboratory, Department of Physics, 
University of Oxford, Oxford, 0X1 3PU, United Kingdom 
(February 1, 2008) 



O 
O 

» 

o 
O 

On 
(N 



> 

in 
o 



o 

c3 



c 

o 
o 



We derive quantum evolution equations for the dynamics 
of dilute condensed Bose gases. The approach contains, at dif- 
ferent orders of approximation, for cases close to equilibrium, 
the Gross Pitaevskii equation and the first order Hartree Fock 
Bogoliubov theory. The proposed approach is also suited for 
the description of the dynamics of condensed gases which are 
far away from equilibrium. As an example the scattering of 
two Bose condensates is discussed. 

PACS numbers: 03.75.Fi, 34.50.-s 



I. INTRODUCTION 



The description of the evolution of condensates far 
from equilibrium is now a matter of considerable im- 
portance to the field of matter wave physics. This in- 
cludes colliding jjj and collapsing [|||] condensates where 
components of the overall condensed system are moving 
rapidly with respect to one another. Numerous quanti- 
tative theoretical studies of partially condensed gases to 
date have relied mostly on the time dependent Gross- 
Pitaevskii approach and extended kinetic theories (see, 
e.g., that are appropriate to fluctuations and ex- 

citations around thermal equilibrium situations. These 
general approaches all depend on the contact potential 
approximation of the binary interactions which restricts 
two body collisions to zero momenta. New methods are, 
therefore, required to describe systems outside this re- 
stricted domain. 

In this paper we shall describe a technique that we 
believe can provide a comprehensive approach to the de- 
scription of evolving condensates. We use the method 
of non-commutative cumulants to derive a set of cou- 
pled equations that should provide an accurate descrip- 
tion in regimes where previous approaches fail. We shall 
explain the formal basis of the method before applying 
it to the partially condensed gas. We shall then focus 
on the specific case of a colliding pair of condensates 
to show the utility of the method. The general coupled 
equations should also be useful for the description of four- 
wave mixing experiments, for collapsing condensates and 
other cases where bulk relative motion of condensates is 
observed. In the next section we shall define the Hamilto- 
nian for the systems we wish to study and find the formal 
equations of motion for the correlation functions. 



II. QUANTUM EVOLUTION EQUATIONS FOR 
CORRELATION FUNCTIONS 

A many body system of identical boson atoms with a 
pair interaction V^r) is described by the Hamiltonian 



H = j Gb^t(x)H 1B (x)^(x) 



<ixi 



dx 2 V' t (xi)-!/' t (x2)V'(x 1 - x 2 )^(x 2 )V'(xi) 



(1) 



Here Hib(x) is the one body Hamiltonian containing 
the kinetic energy —ti 2 A x /2m and the trapping potential 
Vtrap(x). The field operators satisfy the boson commu- 
tation relations 



[V(x 1 ),^(x 2 )] = 0, 
?A(xi),'0 t (x2)] =<5(xi-x 2 ) 



(2) 



Throughout this article the state of the many body 
system at time t is described by a statistical operator 
p(t) and the expectation value of an arbitrary operator 
A is denoted by 



(A) t = Tr(p(t)A). 



(3) 



The Schrodinger equation then yields the dynamic equa- 
tion 



ih^{A) t = {[A,H]) t . 



(4) 



In the following, expectation values in Eq. (||) are de- 
noted as correlation functions if A is a product of cre- 
ation or annihilation operators. Since every expectation 
value can be expanded in linear combinations of normal 
ordered correlation functions their dynamics determines 
all physical properties of the many body system. Accord- 
ing to Eq. ([|) the system of dynamic equations for the 
normal ordered correlation functions up to the order of 
n reads 



at 

d 

^^V0(x2)^(xi))t = 
iftj^t(x 2 Mxi)) t 



;[^(x 2 )v(x 1 ),if]) t , 

([V- t (x 2 )V(x 1 ),i/]> 4 , 



d 

ih—(lp'(Xr, 



•V(x 1 )) 4 = <[^ t (x„)-->(x 1 ),iJ]) t . 



(5) 
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The system Eq. (5) is incomplete since the commu- 
tator of a product of n field operators with H contains 
products of n + 2 field operators. The main problem 
thus consists in truncating Eq. (5) appropriately. An ap- 
proximate method to close the kinetic equations at an 
arbitrary order has been proposed by Fricke (7| in or- 
der to describe the dynamics of many body systems on 
short time scales. In the latter approach Eq. (5) is trans- 
formed into an equivalent system of differential equations 
for what are termed non-commutative cumulants. This 
article contains the extension of Fricke's work to the de- 
scription of the dynamics of dilute Bose gases and their 
equilibrium properties. 



III. DYNAMIC EQUATIONS FOR NON 
COMMUTATIVE CUMULANTS 

Given a set of boson creation and annihilation opera- 
tors Bi, B%, . . . the cumulant expansion of their respec- 
tive correlation functions is defined recursively by M 

(Bi) = (5i) c , 

{B 1 B 2 ) = {B 1 B 2 ) c + {B 1 ) c {B 2 y i 
{B X B 2 B 5 ) = {B l B 2 B 3 Y+{B l ) c {B 2 B z ) c 

+ (B 2 > C (B 1 B 3 > C + (B 3 ) C (S 1 S 2 ) C + {B l ) c {B 2 Y{B 5 ) c , 

'■ (6) 

The expressions {By. . ■■B n ) c define the non commutative 
cumulants we shall use. The cumulants can be obtained 
directly in terms of correlation functions through 



(5i •••£?„) c 

= _d a_ log , 

dx\ dx n 



.x-iB-i 



) l 211=0,. 



(7) 



In contrast to their respective correlation functions the 
cumulants should decrease with increasing order, as long 
as the system is not too far away from the interaction free 
equilibrium. More precisely, for an ideal Bose gas in ther- 
mal equilibrium all cumulants containing more than two 
field operators vanish. In fact, if the state of an ideal gas 
is given by the grand canonical density matrix, accord- 
ing to Wick's theorem of statistical mechanics || , every 
number conserving normal ordered correlation function 
{B\ ■ ■ ■ B n ) splits into a sum of all possible products of 
contractions (BiBj) which preserve the initial order of 
the operators. Since the expectation value of a single field 
operator vanishes, the contractions equal their respective 
second order cumulants. Inserting the Wick decomposi- 
tions of the correlation functions into the left hand sides 
of Eqs. (6) shows successively that only the number con- 
serving cumulants of the second order can be different 
from zero. In this sense the higher order cumulants are 
a measure of the deviations of the physical system from 
the interaction free equilibrium. 



Among other general properties it is worth mentioning 
that in the absence of a pair interaction the cumulants 
satisfy the same dynamic equations as their respective 
correlation functions. In addition, all field operators in 
(Bi ■ ■ ■ B n ) c commute as soon as the order n exceeds two, 
i.e. 



BiB 



i-Di+l 



Bi+\Bi 



(8) 



According to Eq. (7) the cumulants of the order of n > 1 
depend only on the "deviation operators" Bi — {Bi), i.e. 



{B ± ■ ■ ■ B n ) c = ((Bi - (B t )) ■ ■ ■ {B n - {B n ))Y 



(9) 



Hence, the cumulant expansion Eq. (6) is also useful for 
the description of Bose condensed systems with non num- 
ber conserving correlation functions. 

Due to their particular roles in the description of con- 
densed dilute Bose gases throughout this article, the first 
and second order normal ordered cumulants are denoted 
as follows: 



*(x,f) = <V(x)>£, 
$(xi,x 2 ,<) = (^(x 2 )V(xi))£, 

r(x 1 ,x 2! t) = (^t( X2 )^( Xl ))c, 



(10) 



The first order cumulant ^(x, t) is usually referred to 
as the condensate wave function and $(xi,x 2 ,t) as the 
pair function. r(xi,x 2 ,i) is interpreted as the one body 
density matrix of the non condensed fraction. The total 
one body correlation function then reads 

p( 1 )(x 1 ,x 2 ,i) = (^ t (x 2 )^(x 1 )) t 

= r(xi,x 2 ,t) + *(xi,t)**(x a ,t). 

(11) 

Hence, the density of the gas at the position x and time 
t is given by p^'{x,x,t). 

The dynamic equations for the normal ordered cumu- 
lants can be derived from Eqs. (5) by commuting succes- 
sively all creation operators in the commutators on the 
right hand sides to the left and expanding the correla- 
tion functions on both sides according to Eq. (6). Hence, 
the exact dynamics of the condensate wave function is 
described by 

i/h|*(x,f) = fr 1B ( x )*(x,t) 

at 

+ J dyy(x-y)[(V t (y)^(y)V(x)) t c + *(y,t)r(x,y,t) 

+ *(x,t)r(y,y,t)] 

+ / dyV{x - y)**(y, t) [$(x, y, t) + $(x, t)tt(y, t)] . 
J (12) 

In the same way the exact equations of motion for the 
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pair function 

d 

lh—$(x 1 ,X 2 ,t) = iJ 2 B(xi,X 2 )$(xi,X2,t) 

at 

+ y(xi-x 2 )*(xi,t)*(x 2 ,t) 

dyy( Xl - y)[(^ t (y)V(y)V'(x 2 )V(xi))? 

+ *(y,t)(^ t (y)V(x 2 )V(x 1 ))^ 
+ *(x 1 ,t)(^( y ) V ,( y ) V ,( X2 ))- 
+ T(y, y, t)$( Xl , x 2 , t) + r(xi, y, i)$(x 2 , y, i)] 

+ /" dyy(x 1 -y)**(y,t)[(^(y)^(x 2 )V(xi)) t c 

+ *(y, t)$(xi , x 2 , i) + *(xi, t)*(x 2 ,y, i)] 

+ /" 4yV(xi-y)r(x 2 ,y,t) 

x [$(xi,y,i)+tf(xi,i)*(y,i)] 
+ { Xl ^x 2 } (13) 
and the density matrix of the non condensed fraction 



ih—T(xi,x 2 ,t) = [If ib (xi 



-ffiB(x 2 )]r(xi,x 2 .i) 
dyV( Xl - y)[(^( X2 )V t (y)^(y)^(x 1 )>? 
+ *(y,t)<V t (x2)V t (y)V'(xi))? 

+ *(x 1 ,t)(v>t(x2)V t (y)V(y)) i c 

+ T(y, x 2) t)r(xi, y, i) + r(xi, x 2) t)r(y, y, t)] 

dy^(x 1 -y)**(y,i)[(V t (x2)^(y)^(x 1 )>^ 
+ *(xi,t)r(y, x 2 , t) + (y, t)r(xi, x 2 , t)] 
+ y dyF(x 1 -y)$*(y,x 2 ,t) 

x [$(xi,y,t) + *(xi,t)*(y,t)] 

- {x x ~ x 2 }* (14) 

are obtained. Here {xi <-> X2} denotes the exchange of 
the coordinates xi and X2 in the bracket and 

i?2B(xi,x 2 ) = 7?ib(xi) + J?ib(x 2 ) + V(xi - x 2 ) 

(15) 

denotes the two body Hamiltonian. 

As proposed by Fricke , a first approach to close the 
dynamic equations for cumulants could consist in keeping 
all cumulants up to a certain order n. If, for example, n = 
1 is chosen Eq. (12) yields the Gross Pitaevskii equation 
in the Born approximation, 



ffis(x)- 
x #(x,t). 



dyV(x-y)\V(y,t)f 



Neglecting all cumulants containing three or four field op- 
erators in Eqs. (12), (13) and (14) yields a system of time 
dependent first order Hartree Fock Bogoliubov equations 
which constitute the second order approximation. An 
alternative derivation of these second order equations, 
in the contact potential approximation, which includes 
applications to the Timmermans model of molecular for- 
mation in condensates || can be found in Ref. [lCf] . An 
independent approach which contains also the cumulants 
up to the third order has been proposed by Proukakis and 
one of the authors [jll| . 

These methods neglect multiple scattering in the equa- 
tions for the last two relevant cumulants and can be in- 
terpreted as an expansion in terms of V(t — t )/h, where 
t — to is the time passed since the initial time. For this 
reason they are applicable mainly for the description of 
short time dynamics fit] . The multiple scattering contri- 
butions are contained in the dynamic equations for the 
cumulants of the order of n + 1, n + 2, etc.. 

A systematic extension of the previous methods, which 
accounts for multiple scattering in all equations of mo- 
tion for the normal ordered cumulants up to the order 
of n, can be achieved by including the free dynamics of 
the normal ordered cumulants of the order of n + 1 and 
n + 2. Here, free dynamics means neglecting all products 
of normal ordered cumulants containing n + 3 or n + 4 
field operators in the equations of motion for the nor- 
mal ordered cumulants of the order of n + 1 and n + 2. 



As will be shown in Subsection IVA, including multiple 



scattering in all equations of motion at a given order of 
approximation significantly extends their range of valid- 
ity. 



IV. FIRST ORDER APPROACH 

In this section the extended cumulant method will be il- 
lustrated for the lowest order n = 1 approximation. 
As in the case of the related Gross Pitaevskii approach, 
the resulting lowest order dynamic equation is applicable 
only as long as the influence of the non condensed frac- 
tion on the condensate can be neglected. This, in turn, 
restricts the initial states of the gas we can treat with 
this approximation. 

Neglecting all products of cumulants containing four 
or five field operators in the equations of motion for the 
cumulants of the order of two and three, the dynamic 
equation for the condensate wave function is given by 
Eq. (12). According to Eqs. (13) and (14) the respective 
equations for $ and T in the same approximation become 

d 

ih— $(xi,X2,f) =i/ 2B (xi,x 2 )$(xi,x 2 ,i) 
at 



V(xi-x 2 )*(xi,t)*(x 2 ,t) 



(17) 



(16) 
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and 



iH— r(xi,x a ,t) = [7J 1B (x 1 )- J ff 1B (x 2 )]r(x 1 ,x 2 ,i). 
dt (18) 

Among all third order normal ordered cumulants, 
only those which appear on the right hand side of 
Eq. (12) have to be considered. We first expand 
the dynamic equation for the correlation function 
(tp^(x3)ip(x2)ip(xx))t according to Eq. (6) and neglect 
products of normal ordered cumulants containing five 
field operators. The equation of motion for the only rel- 
evant normal ordered third order cumulant is then given 
by: 

^(^ t (x 3 )V(x 2 )^(x 1 ))^ 

= [ff 2B (x 1 ,X 2 ) - ff 1B (x3)](V t (x 3 )V(x 2 )V(x 1 )) t C 



7(xi-x 2 )[*(x2,t)r(xi,xa,t) 

+ *(xi,t)r(x a ,xa,t)]. 



(19) 



Equations (17), ( [Tq ) and (19) can be solved formally in 
terms of the time dependent two body Green's function. 
A detailed description of the Green's function methods 
applied in this article is given in Ref. 0|. The retarded 
two body Green's function is related to the respective 
time development operator t/ 2B (i,io) = exp[— i7J 2B (i — 
to)/h] through 



^2B (*' *o) = _ *o)^2b(*)*o)) 



1 



ih 



(20) 

where 9(t — to) is the step function which yields unity for 
t > to and vanishes elsewhere. Equation (J2fj) shows that 

G 2 g^ obeys the equation: 



ih- 



dt 



2B 



Gg\t,t ) = S(t-t ). (21) 



A convenient representation of the formal solutions to 
Eqs. (17), (|l8|) and (19) is achieved by changing from 
position space to the basis defined by the trap states 
(pi. In this representation the field operators ip become 
single mode annihilation operators ai = J*dx<^(x)^>(x), 
and the relevant cumulants are given by = (a,i}%, 

$y(t) = {a jai ) c t and r tf (t) = (a) ai ) c t . Using Eq. © the 
formal solution for the pair function in Eq. (17) can be 
written in the form: 

"MO = ^(i, i|t^2B(i, to)|fci, fc 2 > s ^>fc 1 fe 2 (*o) 
feifea 

+ [ dt 1 .(iJ\G&\t,t 1 )V\k 1 ,h a ) B * hl (t 1 )* ka (t 1 ). 
kik 2 Jt ° (22) 

Equation (|l8] ) yields the density matrix of the non con- 
densed fraction thus: 



In the same way, the formal solution for the third order 
cumulant in Eq. (19) becomes 



(a\aja,i)l = ^ s(i, j\U 2 B(t, t )\ki, k 2 ) s 

k\k 2 

k 1 k 2 Jt ° 

x ^k 1 {ti)Tk 2 k{tiy Ek(t ' tl),h - (24) 

Here to denotes the initial time and \i,j) s is the (not nec- 
essarily normalized) symmetrized product of trap states 
[0i(xi)0j(x2) + 0i(x 2 )<#j(x:i)]/2. According to Eqs. @ 
and (24) in the first order cumulant approach, the den- 
sity matrix of the non condensed fraction (t) and the 
third order cumulant (a k ajai)f will not evolve in time 
as long as they vanish initially. Their contributions can 
therefore be neglected. 

Equations (22), (§|) and (24) inserted into Eq. (12) 
yield the closed nonlinear Schrodinger equation for the 
condensate wave function 

/CO 
dt ls {i,k 3 \T^+\t,t 1 )\k 1 ,k 2 ) s 



(25) 



Here, T 2B denotes the time dependent and retarded two 

(26) 



body transition matrix 

T^\t,t ) = VS(t - to) + VG 2 l\t,t Q )V. 



T ij {t) = e- i ^- E ^ t - t ^' h T ij {t Q ). 



(23) 



The initial pair function $fc 1 fc 2 (^o), which accounts for 
collisions between condensate atoms that occur before 
t , is neglected in Eq. (25). If the condensate is prepared 
as an ideal gas ground state &kik 2 (to) vanishes exactly 
Jl3| . Recent experiments using magnetic fields to manip- 
ulate the pair interaction fl4|| have shown how to prepare 
such an initial state of a Bose condensed gas. In general, 
the state of the gas will include these initial correlations 
produced by collisions before the initial time we are con- 
sidering. 

It is worth noting that according to Eqs. (^6|) and ( |20| ) 
in the time integral on the right hand side of Eq. (25) the 
condensate wave function is evaluated only at the present 
time t\ = t or in the past. This so called non-Markovian 
property of the nonlinear Schrodinger equation results 
from the influence of the higher order cumulants on the 
condensate, which modifies its dynamics. 



A. Equilibrium properties 

To describe many present day experiments we can use the 
Markov approximation in Eq. (25). In order to obtain 
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the Markov limit, the condensate is assumed to be in an 
ideal gas initial state and the condensate wave function 
is transformed to the interaction picture 

= ^\{t)e- iEit ' h . (27) 

For a dilute gas the non Markovian interaction term on 
the right hand side of Eq. (25) can be considered as a 
small perturbation of the first, i.e. ideal gas, contribution. 
This in turn leads to a weak time dependence of ^\{t). 
In the interaction picture representation Eq. (25) reads 



8 f°° 
ift— = ]T / dt ls (i,k 3 \T : 

m k lk2k3 Jt ° 



£\t,t 1 )\k 1 ,k 2 ) 



(28) 

According to Eq. ( ]26| ) the time dependent transition ma- 
trix (t,ti) is sharply peaked at ti = t with a width 
determined by the two body collisional duration. As 
^\(t) is slowly varying on this time scale, all interac- 
tion picture condensate wave functions can be evaluated 
at t% = t in Eq. (28). The remaining time integral in 
Eq. (28) then simply contains a Fourier transform of the 
time dependent transition matrix Eq. ( p6[ ) and plays the 
role of a coupling function 

/t — 1 
dT s (i,k 3 \T^'(t, t - T)|fci,fc 2 ) s 
-OO 

x e i(E kl +E k2 )r/h_ (29) 

Here, T 2B becomes approximately independent of time as 
soon as t ~ to exceeds the two body collisional duration. 
According to Eq. (|2^) the time dependent transition ma- 
trix T!^" ' is related to the usual two body T matrix [ |l2] ] 
through 

T^\t, t - r) = -L / dEe- lET ' h T 2B (E + iO), 

2ttH J (30) 

which still contains the trapping potential. Here "iO" 
denotes an imaginary energy "ze" , where the positive pa- 
rameter e is taken to zero after performing the energy 
integral. In many present day experiments the trapping 
potential is nearly constant on the scale of the spatial 
range of the two body interaction potential, so the two 
body T matrix in Eq. ( |30| ) can be replaced by its free 
space counterpart T^ e ■ The time limit t — to — > 00 in 
Eq. (29) then yields 

T 2B {i,k 3 ;ki,k 2 ) = s (i,k 3 \T^ c (E kl + E k2 +i0)\k 1 ,k 2 ) s . 

(31) 

As a consequence, the condensate wave function ^i(t) 
Eq. (28) obeys the following form of the Gross Pitaevskii 
equation: 



at 



+ }^ T 2B (i,k 3 ;k l7 k 2 )^ kl (t)^ k2 (t)^* k3 (t). 

k 1 k 2 k 3 



The coupling function in Eq. (|3l|), however, exhibits an 
imaginary part which describes collisional losses of ener- 
getically excited condensate atoms. The physical signifi- 
canc e of t his imaginary part will be explained in Subsec- 

ivb| . 



tion 

Close to equilibrium, the energy of occupied modes can 
be considered to be so small that the coupling function 
is purely real to a good approximation. The two body 
T matrix can then be evaluated by means of the contact 
potential approximation 

T 2B e = J dp 1 dp 2 dp 3 dp A S(p 4 + pa - Pi - P2) 



4nh 2 

x |P4>|P3> a (pi\(p 2 \ 

m 



(33) 



Here |pj) denotes a one body plane wave momentum 
state and m the atomic mass while ao is the two body s 
wave scattering length. In the contact potential approx- 
imation Eq. (32) can be shown to be formally equivalent 
to the Gross Pitaevskii equation in the Born approxi- 
mation (see Eq. (16)) with the two body potential V(r) 
replaced by 47rfi, 2 ao<5(r)/m. Hence, in the position repre- 
sentation, the equilibrium limit of the Markovian nonlin- 
ear Schrodinger equation in Eq. (32) has the form: 



Hib(x) 



4irh 2 



a |*(x,t)| s 



*(x,t). 
(34) 



(32) 



This Gross Pitaevskii equation is transformed to its time 
independent form by means of the ansatz ^(x, t) = 
\p(x) exp{—ifMt/h), where the energy parameter /1 is usu- 
ally termed the chemical potential. The Markov approx- 
imation of the first order cumulant approach thus leads 
to the correct Gross Pitaevskii equation, which describes 
virtually all equilibrium phenomena in present day ex- 
periments with Bose condensates near zero temperature. 



B. Scattering of two condensates 

In this subsection, the cumulant approach will be ap- 
plied to describe the scattering of two condensates: a 
clear and important non-equilibrium case. In a typical 
experiment the different condensates are generated from 
a trapped parent condensate, at nearly zero temperature, 
by means of an optical standing wave pulse. The trap in 
some experiments can be switched off directly after the 
application of the light pulse at time to. The photon 
momentum p p h can be transferred to the atoms and is 
assumed to be much larger than the momentum spread 
in the parent condensate and also the momentum mc as- 
sociated with the speed of sound c. This implies that 
the condensate wave function prepared in this way con- 
tains three well separated momentum components with 
the respective central atomic momenta p = ±2p p h and 
p = 0. For simplicity, the following analysis is restricted 
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to the presence of two distinct components. The conden- 
sate wave function can then be represented in the center 
of mass frame as [[151 



*(x,to) = * + (x,io)e lPph - x/n + *-(x,i )e- 4p p h - x/ri 

2 



x e 

Here the wave functions 

* + (x,t ) = A + *(x), 

and 

tf_(x,to) = A_*(x), 



(35) 
(36) 
(37) 



are multiples of the stationary wave function of the par- 
ent condensate Vt'(x). The latter is a solution of the 
time independent Gross Pitaevskii equation containing 
the trapping potential. The amplitudes A + and A_ are 
normalized as \A + \ 2 + |A_| 2 = 1. This ensures that 
/ <£k|\P(x, io)| 2 equals the number of atoms in the parent 
condensate, as long as the overlap in momentum space 
of the two components on the right hand side of Eq. (35) 
is negligible. 

The following analysis of the scattering process is based 
on our non-Markovian nonlinear Schrodinger equation 
given by Eq. (25). In position space this equation has 
the following form: 

ifi|U(x,t) =#i B (x)*(x,t) 

/'CO 

dyidy 2 dy 3 / *i{x,y 3 |T 2 ( B ) (t ) ti)|y 1 ,y 2 ) 

J t 

x*(yi,ti)*(y 2 ,ti)**(y 3j t)- (38) 

The initial condition Eq. (35) implies that the represen- 
tation 



*(x,i) = [tf + (x,t)e ip <*- x /' t + *_(x,t)e- ip ^- x / fi 



x e 



-i^h. t /H 



(39) 



provides a suitable interaction picture, where \& + (x, t) 
and \P_(x, i) are assumed to be slowly varying in space 
and time. This interaction picture representation was 
proposed by Band et al. in Ref. Jlj| and is equivalent 
to the slowly varying envelope approximation SVEA of 
nonlinear optics flqj . Equation (39) inserted into Eq. (38) 



yields the following equation of motion: 

h 2 



P P h 

rn 



Vx | + ^A, 



*+(x,i) 



+ {p P h <-> -Pph} 

/f°° p2 i 
dyidy 2 dy 3 / dtie^*-^ 

x (x,y 3 |T 2 ( + ) (i,t 1 )|y 1 ,y 2 > 

x r e l Pph-yi/ R ^ + ( yi ,< 1 ) + e -*pp"-yi/'^_( yi ,t 1 ) 

x [e 4p ^' y2/R * + (y 2 ,<i) +e- lp p h - y2 / fi *_(y 2 ,t 1 ) 
x \e~ lp ^- y3 / h ^* + (y 3 ,t) + e tp ^ y3 / h y*_{y 3 ,t) 

Here, (t, t\) is the time dependent transition oper- 
ator Eq. (|26J) in free space and {p p h *-* ~P P h} means 
the exchange of p p h and — p p h and of ^ + and in- 
side the bracket. The interaction picture wave functions 
^±(yi ; ti) an d v P±(y2, t±) are slowly varying on the time 
scale of the two body collisional duration and can be eval- 
uated at t% = t on the right hand side of Eq. (40). This 
leads to the Markov approximation. In virtue of Eq. j30| ) 
the remaining time integral yields 



(40) 



(41) 



as soon as £ — to exceeds the two body collisional duration. 

To perform the spatial integrations on the right hand 
side of Eq. (40) the position representation of the T ma- 
trix in Eq. pi] ) is 
and becomes 



transformed to the momentum space 



(x,y 3 |T 2 B(Pp h /TO + «0)|yi,y 2 ) 

^— g J dp 1 dp 2 dp 3 dp i S{p 1 + p 2 

'Pph (Pi +P2) 2 



(2-Khy 

P3 P4 



21 



4m 




x e »(P4-x+p3-y 3 )/?i e -i(pi-yi+p2-y2)/?! 



where T 2 b is the T matrix which contains only the rela- 
tive motion of pairs of atoms. The momentum conserving 
5 function stems from the translational invariance of T 2 b 
in the center of mass coordinates. For the experimen- 
tally relevant photon momenta p p h, the T matrix on the 
right hand side of Eq. (42) depends significantly on the 
related kinetic energy p 2 h /m and the contact potential 
approximation Eq. (33) cannot be used. In the Markov 



G 



approximation Eq. (40) then reads 



the following equation of motion: 



, , 9 p ph 

at m 



2m 



-A, 



*+(x,t) 



{Pph -Pph} 

==3 / dpidp 2 dp 3 rfp4*(Pl + P2 P3 - P4) 

Pi P2 




?2B — j h ?0 



x \ [*+(pi - Pph, t)tf + (pa - Pph, t)*+(p3 - Pph, t) 
*+(Pi - Pph, t)*-(pa + Pph, t)*+(pa - Pph, t) 

*-(Pl +Pph,t)*+(P2 -Pph,*)*+(P3-Pph,t) 
*+(Pl - Pph, *)*+(P2 - Pph, t)**(p3 + Pph, t)] 



+ [p P h--Pph]}e ip4 - x/ri , 
where ^±(p, t) denotes the Fourier transform 

*±(P,*) = ^^3 / dxe-W^^t). 



(43) 



(44) 



The wave functions ^±(x,t) are slowly varying in 
space, and so their Fourier transform v E'±(p, t) is sharply 
peaked at zero momentum (p = 0) with a width de- 
termined by the initial momentum spread of the par- 
ent condensate. According to Eq. (p0|), the two body 
T matrix T 2 b(E + iO) is slowly varying on the scale of 
the energies corresponding to the inverse two body col- 
lisional duration. Its plane wave basis matrix elements 
(Pout|32s|Pin) are slowly varying in p in and p ou t, on the 
scale of the momentum associated with the inverse range 
of the two body interaction potential V(r). This can be 
seen through the relation: 



(rout|72 B |rin) =V(r in )S(r in - r out ) 

+ V r (r ou t)(r out |G 2 B|rin)V r (r in ), 



(45) 



where G 2 b denotes the two body energy-dependent 
Green's function for the relative motion only. As a con- 
sequence, the plane wave matrix element of the T matrix 
T 2 b on the right hand side of Eq. (43) can be evaluated 
at the respective momenta determined by the products of 
momentum space condensate wave functions and by the 
momentum conserving 5 function. Performing the mo- 
mentum space integrals then transforms the condensate 
wave functions ty± back to position space and leads to 



at m 



2m 



+ {Pph <-> -Pph} 



a ip p h-x/R 



ink 2 



o |*+(x,t)| 2 *+(x,t) 



[/(k ph , k ph ) + /(k ph , -k ph )] |*_(x, t)\ 2 *+(x, t) 



+ e 3lp ^ h <-2p ph \f 2B (iO) |0}¥* (x, t)9% (x, t) 



{p P h <-> -P P h} • 



(46) 



Here, k p h = p p h/^ denotes the wave vector associated 
with the photon momentum p p h and 



/(k out ,k in ) = - 2-K 2 mh 

x (hk out \f 2B (h 2 k 2 /m + i0)\hk h 



(47) 



is the scattering amplitude, where the wave vectors k; n 
and k out are on the energy shell |k out | = |k in | = k. Equa- 
tion (46) contains two well separated momentum compo- 
nents on its left hand side and four on the right hand side. 
Keeping only the phase matched terms yields a system 
of two coupled nonlinear Schrodinger equations. These 
are given by 



■ I d fc k Ph 

m 



dt 
Airh 2 

m 
4irh 2 

m 



ooi*+(x,t)r 



^ A 

2m 



*+(x,t) 



[/(kph, kph) + /(kph, -kph)] 



x |vj/_(x,i)| 2 *+(*,*) 



(48) 



for the wave function and the analogous relation for 
1 J r _ where k p h and — k p h as well as ^+ and are ex- 
changed. 

As a result of the Markov approximation the T ma- 
trix elements in Eq. (48) are evaluated on the two-body 
energy shell and can be represented in terms of the scat- 
tering amplitude Eq. (47). The experimentally relevant 
photon momenta are sufficiently small for the scattering 
amplitudes on the right hand side of Eq. (48) to be in 
the isotropic limit. According to effective range theory 
@ /(k ou t, ki n ) in Eq. (47) can then be expanded in the 
form: 

/(kout, ki„) = -ao + iao(kao) + 0(k 2 ). (49) 



Equation (|49j) inserted into Eq. (48) recovers the result 
of Band et al. in Ref. |15| ). In the latter article, the two 
body energy conservation in the collisions was assumed 
for physical reasons and the imaginary part on the right 
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hand side of Eq. (|49|) included, on an heuristic basis in 
the Gross Pitaevskii equation, in the form of Eq. (34). 



The imaginary correction term includes the loss of those 
atoms from the condensates that are elastically scattered 
out of the forward direction. As reported in Ref. [15|, 
these losses have been observed experimentally in four 
wave mixing experiments with Bose condensates jlj and 
exhibit an excellent agreement with the predictions of the 
appropriate generalization of Eq. (48). 
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V. CONCLUSIONS 



We have shown in this paper how we may systematically 
obtain a theoretical description of the dynamics of dilute 
Bose gases which are far from their thermal equilibrium. 
The underlying cumulant approach is based on truncat- 
ing the infinite hierarchy of quantum evolution equations 
for correlation functions of boson field operators in ac- 
cordance with Wick's theorem of statistical mechanics. 
We extend this approach by giving a general scheme for 
including multiple scattering in all dynamic equations. 
This extended approach describes the inter-particle colli- 
sions in terms of scattering amplitudes rather than bare 
potentials. 

The lowest order approximation consists in a non 
Markovian nonlinear Schrodinger equation. When ap- 
plied to initial conditions close to the thermal equilibrium 
at zero temperature our lowest order approach recovers 
the time dependent Gross Pitaevskii equation. Applied 
to the scattering of two condensates at a sharply defined 
relative velocity, which is high in comparison with the in- 
ternal velocity spread in each condensate, we have shown 
that the non Markovian nonlinear Schrodinger equation 
goes beyond the Gross Pitaevskii approach and describes 
the experimentally observed collisional loss of condensate 
atoms. Hence, including multiple scattering not sim- 
ply justifies the common contact potential approximation 
but accounts for important new physical phenomena. 

The second order basic cumulant approach in the con- 
tact potential approximation was shown to be equivalent 
to the first order time dependent Hartree-Fock Bogoli- 
ubov theory derived in Ref. j[o|. A detailed description 
of the physical significance of the corresponding extended 
approach, which includes multiple scattering, goes be- 
yond the scope of this article and will appear elsewhere 

Our first order approach should also be applicable 
when the relative bulk motion is still rapid but its veloc- 
ity is not sharply defined. This situation occurs during 
the collapse of a Bose condensed gas with a negative s 
wave scattering length. Very recently, the dynamics of a 
collapse and the related losses of condensate atoms have 
been determined experimentally in Ref. |3|. 
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